/******************************************************************************
Author: Akshay Dixit
Date Created: May 23, 2018
Date last modified: June 11, 2020
File Name: 4_in.do
Project: T4D Indonesia 
Purpose: Analysis of intermediate outcomes (Household data)
******************************************************************************/

/* LIST OF INTERMEDIATE OUTCOMES & THE ASSOCIATED VARIABLES (from the household dataset)
	
	Note that the following codes apply:
		V. Other, specify _________
		X. REFUSED
		Y. DON'T KNOW
		
		95. Other, specify _________
		97. REFUSED
		98. DON'T KNOW
		99. MISSING

Category 1. Increased awareness, knowledge and improved community  attitudes

	Outcome 1:	Number of questions on attitudes towards healthcare practices answered appropriately by the respondent. 
				sp01 sp02 sp03
	Outcome 2:	Number of healthcare details that the respondent was able to correctly list.
				pf08 pf10 pf41_d 
	Outcome 3:	Number of questions on knowledge of health levers that the respondent is able to correctly answer.
				sp06 sp06_a sp06_b sp06_c sp07 sp08 sp09_a sp09_b
	Outcome 4:	Number of questions on partner participation/decision making that the respondent is able to answer appropriately.
				sp04 sp05 sp10 pf16a
				
Category 2. Improved facility access (transportation, new facility, longer facility hours, outreach services)
				
	Outcome 2:	Community organized transportation
				pf30 td04_dumL1
	Outcome 3:	Fix road [see .do file titled "section_td"]
				td01H 
	Outcome 4:	Travel time
				pf31_x pf31_a pf31_b
	Outcome 5:	Travel cost
				pf32 pf33

Category 5. Increased ability to pay

	Outcome 1:	Address cost of services - Whether the respondent reported cost of care as a barrier to utilization of care.
				pf24
	Outcome 2:	Address cost of services - Total fees (if any) for the delivery reported by the respondent
				pf27_x pf27
	Outcome 3:	Address cost of services - Whether the respondent (or spouse) has an insurance/health protection program enrolment
				pf27a_x pf27a pf28 pf28a
	Outcome 4:	Raise community funds - Whether the respondent was part of a woman's savings group to help save for costs associated with MNH care.
				pk17
	Outcome 5:	Raise community funds - Whether the respondent reported using a community fund to help pay for costs associated with MNH care.
				pk18 pk19

Category 9. Improved facility cleanliness
	
	Outcome 1:	Level of cleanliness reported by the respondent
				pk12A pk12B

Category 10. Improved attitude, effort, trust of the provider		
	
	Outcome 1:	The level of satisfaction reported by the respondent with the quality of MNH services
				pk07* pk08* pk09* pk10* pk11* pk13* pk14* pk15*
	Outcome 2:	If the respondent reported non-dignified care. 
				pf34 pf35 
	Outcome 3:	If the respondent reported provider neglect.
				pf36 pf37
	Outcome 4:	If the respondent reported provider abuse. 
				pf38
*/

clear all
set more off

u "$data/household_data.dta", clear

********************************************************************************

* CATEGORY 1. Increased awareness, knowledge and improved community  attitudes

***Preparing Outcome 1***
desc sp01 sp02 sp03
foreach var in sp01 sp02 sp03 {
	replace `var' = . if `var' == 97
	g `var'_binary = (`var' == 3 | `var' == 4)
	replace `var'_binary = . if `var' == .
	tab `var'
	tab `var'_binary
}

egen knowledge_practices = rowtotal(sp01_binary sp02_binary sp03_binary), missing
replace knowledge_practices = . if sp01_binary == . | sp02_binary == . | sp03_binary == .
	//Only include complete cases in the analysis
tab knowledge_practices

***Preparing Outcome 2***
desc pf08 pf10 pf41_d 
foreach var in pf08 pf10 pf41_d {
	tab `var'
}

*PF08
foreach value in A1 B1 C1 D1 E1 F1 G1 H1 I1 J1 K1 L1 M1 N1 O1 P1 Q1 R1 S1 T1 U1 A2 B2 C2 D2 E2 F2 G2 H2 I2 J2 K2 L2 M2 N2 O2 P2 Q2 R2 S2 {
	g pf08_`value' = strpos(pf08, "`value'") > 0 
}
egen knowledge_problems = rowtotal(pf08_A1-pf08_S2), missing

local problems pf08_A1 pf08_B1 pf08_C1 pf08_D1 pf08_E1 pf08_F1 pf08_G1 pf08_H1 pf08_I1 pf08_J1
foreach var of local problems {
	tab `var'
}

local problems pf08_K1 pf08_L1 pf08_M1 pf08_N1 pf08_O1 pf08_P1 pf08_Q1 pf08_R1 pf08_S1 pf08_T1 pf08_U1 ///
pf08_A2 pf08_B2 pf08_C2 pf08_D2 pf08_E2 pf08_F2 pf08_G2 pf08_H2 pf08_I2 pf08_J2 pf08_K2 pf08_L2 pf08_M2 pf08_N2 ///
pf08_O2 pf08_P2 pf08_Q2 pf08_R2 pf08_S2
foreach var of local problems {
	tab `var'
}
 
*drop pf08_K1-pf08_S2

*PF10
foreach value in A B C D E F G {
	g pf10_`value' = strpos(pf10, "`value'") > 0 
}
egen knowledge_preparedness = rowtotal(pf10_A-pf10_G), missing
drop pf10_A-pf10_G 

*PF41_d
foreach value in A B C D E F G H I J K L M {
	g pf41_d_`value' = strpos(pf41_d, "`value'") > 0 
}
egen knowledge_complications = rowtotal(pf41_d_A-pf41_d_M), missing

local complications pf41_d_A pf41_d_B pf41_d_C
foreach var of local complications {
	tab `var'
}
 
*drop pf41_d_D-pf41_d_M

***Preparing Outcome 3***
desc sp06 sp06_a sp06_b sp06_c sp07 sp08 sp09_a sp09_b

g anc_timing = (sp06 == 1 | (sp06_a <= 12 & sp06_a != .) | (sp06_b <= 3 & sp06_b != .) | (sp06_c == 1 & sp06_c != .))
replace anc_timing = . if sp06 == 97
tab anc_timing

tab sp07
g anc_visits = (sp07 >= 4)
tab anc_visits

g postnatal_yesno = (sp08 == 1)					
g postnatal_timing = (sp09_a <= 7 | sp09_b <= 1) //All conditional missing values and Don't Knows coded as 0

egen knowledge_levers = rowtotal(anc_timing anc_visits postnatal_yesno postnatal_timing), missing
replace knowledge_levers = . if anc_timing == . | anc_visits == . | postnatal_yesno == . | postnatal_timing == .
	//Only include complete cases in the analysis
tab knowledge_levers

***Preparing Outcome 4***
desc sp04 sp05 sp10 pf16a

tab sp04
g sp04_binary = (sp04 == 1 | sp04 == 2)
replace sp04_binary = . if sp04 == 97

tab sp05
g sp05_binary = (sp05 == 3 | sp05 == 4)
replace sp05_binary = . if sp05 == 97

tab sp10
g sp10_binary = (sp10 == 1)
replace sp10_binary = . if sp10 == .

tab pf16a 
replace pf16a = 0 if pf16a == 3

egen partner_participation = rowtotal(sp04_binary sp05_binary sp10_binary pf16a), missing
replace partner_participation = . if sp04_binary == . | sp05_binary == . | sp10_binary == . | pf16a == .
	//Only include complete cases in the analysis
tab partner_participation

egen partner_participation_2 = rowtotal(sp04_binary sp05_binary pf16a), missing
replace partner_participation_2 = . if sp04_binary == . | sp05_binary == . | pf16a == .
	//Only include complete cases in the analysis


***Output Tables***

*Label variables*
lab var knowledge_practices "Knowledge of ANC and facility birth"
lab var knowledge_problems "Knowledge of pregnancy problems"
lab var knowledge_preparedness "Knowledge of birth preparedness"
lab var knowledge_complications "Knowledge of postnatal complications"
lab var knowledge_levers "Knowledge of T4D levers"
lab var partner_participation "Attitudes on partner participation"

*Produce Excel output
cd "$analysis"

//Create the table shell
putexcel set "intermediate_outcomes.xlsx", sheet("Category 1") modify

putexcel A1=("Outcome") B1=("Treatment Mean") C1=("Control Mean") D1=("Impact") ///
E1=("p-value") F1=("Effect Size") G1=("Sample Size") A9=("Number of Respondents") ///
A10=("Number of villages") B9=(3016) C9=(2985) B10=(100) C10=(100) ///
A12=("Treatment means are regression adjusted") A13=("*** p<0.01, ** p<0.05, * p<0.1")

//Export variable label, regression coefficient, control mean, treatment mean, p-value, sample size
//Effect sizes are exported separately (see below)
local outcomes knowledge_practices knowledge_problems knowledge_preparedness knowledge_complications knowledge_levers partner_participation
local row = 2
foreach var of local outcomes {
	sleep 2000
	local varlabel : var label `var'
	
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	
	local p = 2*ttail(e(df_r), abs(_b[treatment]/_se[treatment]))
	local stars
		if `p' < 0.10 local stars = "*" 
		if `p' < 0.05 local stars = "**" 
		if `p' < 0.01 local stars = "***"
	local impact = _b[treatment]
	
	qui sum `var' if treatment == 0
	putexcel A`row' = ("`varlabel'")
	putexcel B`row' = ((_b[treatment]) + (r(mean)))
	putexcel C`row' = (r(mean))
	putexcel D`row' = ("`impact'" + "`stars'")
	putexcel E`row' = (`p')
	putexcel G`row' = (e(N))
	local ++row
}

//Export effect sizes for outcomes that aren't already standardized
local simple_outcomes knowledge_practices knowledge_problems knowledge_preparedness knowledge_complications knowledge_levers partner_participation
local row = 2
foreach var of local simple_outcomes {
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	qui sum `var' if treatment == 0
	putexcel F`row' = (_b[treatment]/r(sd))
	local ++row
}

********************************************************************************

* CATEGORY 2. Improved facility access 

***Preparing Outcome 2***
desc pf30 td04_dumL1
tab pf30
tab td04_dumL1

g travel_on_foot =  strpos(pf30, "F") > 0
replace travel_on_foot = . if pf30 == ""
count if travel_on_foot == 1 & td04_dumL1 == 1		
		//2 instances where respondent travelled both on foot and using community organized transportation

g transport = (travel_on_foot == 0 | td04_dumL1 == 1)
		/* 
		Includes cases where respondents did not travel on foot AT ALL, even in combination
		with other forms of transportation
		
		There are only 2 exceptions to this: i.e. the two cases where respondents
		travelled on foot in combination with community organized transportation
		*/
		
replace transport = . if travel_on_foot == . & td04_dumL1 == .
tab transport

g travel_ambulance = strpos(pf30, "H") > 0 
tab travel_ambulance

***Preparing Outcome 4***
desc pf31_x pf31_a pf31_b
g travel_time = pf31_a
replace travel_time = (pf31_b/60) if pf31_a == . 
tab travel_time

***Preparing Outcome 5***
desc pf32 pf33
tab pf32
tab pf33

	/*PF33 - data on travel cost - is weirdly formatted as text, with numbers
	interspersed with dashes and dots. Fixing this here.
	*/
	count if pf33 != ""	
	forvalues i = 1/9 {
		g pf33_`i' = substr(pf33,-`i',1)
		replace pf33_`i' = "0" if pf33_`i' == "_"
	}
	drop pf33_4 pf33_8
	egen pf33_num = concat(pf33_9 pf33_7 pf33_6 pf33_5 pf33_3 pf33_2 pf33_1)
	destring pf33_num, replace
	drop pf33_1-pf33_9

ren pf33_num travel_cost
replace travel_cost = 0 if pf32 == 3	//If respondent says they didn't spend on travel, travel cost is 0
tab travel_cost	

***Output Tables***

*Label variables*
lab var transport "Traveled to facility by means other than on foot"
lab var travel_ambulance "Traveled to facility by ambulance"
lab var travel_time "Travel time (hours)"
lab var travel_cost "Travel cost (IDR)"


*Produce Excel output

//Create the table shell
putexcel set "intermediate_outcomes.xlsx", sheet("Category 2") modify 

putexcel A1=("Outcome") B1=("Treatment Mean") C1=("Control Mean") D1=("Impact") ///
E1=("p-value") F1=("Effect Size") G1=("Sample Size") A7=("Number of Respondents") ///
A8=("Number of villages") B7=(3016) C7=(2985) B8=(100) C8=(100) ///
A10=("Treatment means are regression adjusted") A11=("*** p<0.01, ** p<0.05, * p<0.1")

//Export variable label, regression coefficient, control mean, treatment mean, p-value, sample size
//Effect sizes are exported separately (see below)
local outcomes transport travel_ambulance travel_time travel_cost
local row = 2
foreach var of local outcomes {
	sleep 2000
	local varlabel : var label `var'
	
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	
	local p = 2*ttail(e(df_r), abs(_b[treatment]/_se[treatment]))
	local stars
		if `p' < 0.10 local stars = "*" 
		if `p' < 0.05 local stars = "**" 
		if `p' < 0.01 local stars = "***"
	local impact = _b[treatment]
	
	qui sum `var' if treatment == 0
	putexcel A`row' = ("`varlabel'")
	putexcel B`row' = ((_b[treatment]) + (r(mean)))
	putexcel C`row' = (r(mean))
	putexcel D`row' = ("`impact'" + "`stars'")
	putexcel E`row' = (`p')
	putexcel G`row' = (e(N))
	local ++row
}

//Export effect sizes for outcomes that aren't already standardized
local outcomes transport travel_ambulance travel_time travel_cost
local row = 2
foreach var of local outcomes {
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	qui sum `var' if treatment == 0
	putexcel F`row' = (_b[treatment]/r(sd))
	local ++row
}


********************************************************************************

* CATEGORY 5. Increased ability to pay

***Preparing Outcome 1***
desc pf24
tab pf24
g cost_barrier =  strpos(pf24, "A") > 0
	//Respondents who reported cost of care as one reason for not delivering at a facility are coded as "1".
	//Respondents who delivered at a facility are coded as "0". 
tab cost_barrier

***Preparing Outcome 2***
desc pf27_x pf27
tab pf27_x
tab pf27

	/* PF27 - data on cost of delivery - is weirdly formatted as text, with numbers
	interspersed with dashes and dots. Fixing this here.
	*/
	count if pf27 != ""	
	forvalues i = 1/10 {
		g pf27_`i' = substr(pf27,-`i',1)
		replace pf27_`i' = "0" if pf27_`i' == "_"
	}
	drop pf27_4 pf27_8
	egen pf27_num = concat(pf27_10 pf27_9 pf27_7 pf27_6 pf27_5 pf27_3 pf27_2 pf27_1)
	destring pf27_num, replace
	drop pf27_1-pf27_10

ren pf27_num delivery_cost
replace delivery_cost = 0 if pf27_x == 3	//If respondent says they didn't spend on delivery, delivery cost is 0
tab delivery_cost
g paid_for_delivery = (delivery_cost != 0)
replace paid_for_delivery = . if delivery_cost == .
	
***Preparing Outcome 3***
desc pf27a_x pf27a pf28 pf28a	
tab pf27a_x
tab pf28
tab pf28a	
g insurance = (pf27a_x == 1 | pf28 == 1 | pf28a == 1)
replace insurance = . if pf27a_x == . & pf28 == . & pf28a == .

***Preparing Outcome 4***
desc pk17
replace pk17 = 0 if pk17 == 3
replace pk17 = . if pk17 == 97
ren pk17 savings_group_member
tab savings_group_member

***Preparing Outcome 5***
desc pk18 pk19
tab pk19
tab pk18
g used_community_fund = (pk19 == 1)
replace used_community_fund = . if pk18 == 97
tab used_community_fund	//If the village does not have a community fund, this variable is 0.

***Output Tables***

*Label variables*
lab var cost_barrier "Cost reported as a barrier to healthcare"
lab var delivery_cost "Cost of delivery"
lab var paid_for_delivery "Paid for delivery (binary)"
lab var insurance "Insurance enrollment"
lab var savings_group_member "Savings group membership for MNH care"
lab var used_community_fund "Used community fund for MNH care"


*Produce Excel output

//Create the table shell
putexcel set "intermediate_outcomes.xlsx", sheet("Category 5") modify 

putexcel A1=("Outcome") B1=("Treatment Mean") C1=("Control Mean") D1=("Impact") ///
E1=("p-value") F1=("Effect Size") G1=("Sample Size") A9=("Number of Respondents") ///
A10=("Number of villages") B9=(3016) C9=(2985) B10=(100) C10=(100) ///
A12=("Treatment means are regression adjusted") A13=("*** p<0.01, ** p<0.05, * p<0.1")

//Export variable label, regression coefficient, control mean, treatment mean, p-value, sample size
//Effect sizes are exported separately (see below)
local outcomes cost_barrier delivery_cost paid_for_delivery insurance savings_group_member used_community_fund
local row = 2
foreach var of local outcomes {
	sleep 2000
	local varlabel : var label `var'
	
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	
	local p = 2*ttail(e(df_r), abs(_b[treatment]/_se[treatment]))
	local stars
		if `p' < 0.10 local stars = "*" 
		if `p' < 0.05 local stars = "**" 
		if `p' < 0.01 local stars = "***"
	local impact = _b[treatment]
	
	qui sum `var' if treatment == 0
	putexcel A`row' = ("`varlabel'")
	putexcel B`row' = ((_b[treatment]) + (r(mean)))
	putexcel C`row' = (r(mean))
	putexcel D`row' = ("`impact'" + "`stars'")
	putexcel E`row' = (`p')
	putexcel G`row' = (e(N))
	local ++row
}

//Export effect sizes for outcomes that aren't already standardized
local outcomes cost_barrier delivery_cost paid_for_delivery insurance savings_group_member used_community_fund
local row = 2
foreach var of local outcomes {
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	qui sum `var' if treatment == 0
	putexcel F`row' = (_b[treatment]/r(sd))
	local ++row
}


********************************************************************************

* CATEGORY 10. Improved attitude, effort, trust of the provider

***Preparing Outcome 1***
desc pk07* pk08* pk09* pk10* pk11* pk12* pk13* pk14* pk15*
ren pk15bA pk15aA
ren pk15bB pk15aB
tab pk01A	//5796 (96.5%) respondents report visiting a facility in the last 12 months
tab pk16A	

	/*
	
	Respondents answered this set of questions for their most recent visit to a facility.
	This most recent visit may or may not have been to a sample puskesmas.
	
	Out of the 5796 respondents who visited a facility in the last 12 months:
		1. For 2095, their most recent visit was to the relevant sample puskesmas.
		2. 1976 respondents did not visit the relevant sample puskesmas at all in the last 12 months.
		3. 1723 did visit a sample puskesmas, but this was not their most recent facility visit.
		
	We analyze here two separate sets of outcomes:
		1. One pertaining to the most recent facility visit (on 5796 observations)
			This includes all respondents who answered this section.
		2. One pertaining only to the sample puskesmas (on 3818 observations - 2095+1723)
			This includes respondents who visited the sample puskesmas in the last 12 months.
			
	*/

local experience pk07A pk07B pk08A pk08B pk08aA pk08aB pk09A pk09B pk10A pk10B pk11A pk11B pk12A pk12B pk13A pk13B pk14A pk14B pk15A pk15B pk15aA pk15aB
foreach var of local experience {
	replace `var' = . if `var' == 98 | `var' == 97
}

g pk07a_A = pk07aA
replace pk07a_A = pk07aA/60 if pk07a_xA == 2

g pk07a_B = pk07aB
replace pk07a_B = pk07aB/60 if pk07a_xB == 2


	*Preparing new variables for analysis pertaining only to the sample puskesmas
	
	forvalues i=1/11 {
		local var: word `i' of "pk07" "pk08" "pk08a" "pk09" "pk10" "pk11" "pk12" "pk13" "pk14" "pk15" "pk15a"
		gen `var'_sample=`var'A if pk16A == 96
		replace `var'_sample=`var'B if pk16A == 1
	}
	
	g pk07a_sample = pk07a_A if pk16A == 96
	replace pk07a_sample = pk07a_B if pk16A == 1
	
***Preparing Outcomes 2,3 and 4***
desc pf34 pf35 pf36 pf37 pf38
local provider_abuse pf34 pf35 pf36 pf37 pf38 
foreach var of local provider_abuse {
	g `var'_binary = strpos(`var', "A") > 0
	replace `var'_binary = . if `var' == ""
	tab `var'_binary
}

g non_dignified_care = (pf34_binary == 1 | pf35_binary == 1)
replace non_dignified_care = . if pf34_binary == . & pf35_binary == .
tab non_dignified_care

g provider_neglect = (pf36_binary == 1 | pf37_binary == 1)
replace provider_neglect = . if pf36_binary == . & pf37_binary == .
tab provider_neglect

ren pf38_binary provider_abuse
tab provider_abuse


***Output Tables***

*Label variables*
lab var pk07A "Recent visit: Satisfaction with waiting time (1-4 scale)"
lab var pk07a_A "Recent visit: Longest waiting time (hours)"
lab var pk08A "Recent visit: Satisfaction with opening hours (1-4 scale)"
lab var pk08aA "Recent visit: Health worker not there (1-4 scale)"
lab var pk09A "Recent visit: Respectful treatment (1-4 scale)"
lab var pk10A "Recent visit: Trust nurses/midwives/other staff (1-4 scale)"
lab var pk11A "Recent visit: Availability of drugs/supplies/equipment (1-4 scale)"
lab var pk13A "Recent visit: Quality of physical facility (1-4 scale)"
lab var pk14A "Recent visit: Overall quality of healthcare (1-4 scale)"
lab var pk15A "Recent visit: If people complained, would staff improve"
lab var pk15aA "Recent visit: Change in the last year (Better care/Worse/Same)"

lab var pk07_sample "Sample Puskesmas: Satisfaction with waiting time (1-4 scale)"
lab var pk07a_sample "Sample Puskesmas: Longest waiting time (hours)"
lab var pk08_sample "Sample Puskesmas: Satisfaction with opening hours (1-4 scale)"
lab var pk08a_sample "Sample Puskesmas: Health worker not there (1-4 scale)"
lab var pk09_sample "Sample Puskesmas: Respectful treatment (1-4 scale)"
lab var pk10_sample "Sample Puskesmas: Trust nurses/midwives/other staff (1-4 scale)"
lab var pk11_sample "Sample Puskesmas: Availability of drugs/supplies/equipment (1-4 scale)"
lab var pk13_sample "Sample Puskesmas: Quality of physical facility (1-4 scale)"
lab var pk14_sample "Sample Puskesmas: Overall quality of healthcare (1-4 scale)"
lab var pk15_sample "Sample Puskesmas: If people complained, would staff improve"
lab var pk15a_sample "Sample Puskesmas: Change in the last year (Better care/Worse/Same)"

lab var non_dignified_care "Non-dignified care" 
lab var provider_neglect "Provider neglect"
lab var provider_abuse "Provider abuse"


*Produce Excel output

//Create the table shell
putexcel set "intermediate_outcomes.xlsx", sheet("Category 10") modify 

putexcel A1=("Outcome") B1=("Treatment Mean") C1=("Control Mean") D1=("Impact") ///
E1=("p-value") F1=("Effect Size") G1=("Sample Size") A28=("Number of Respondents") ///
A29=("Number of villages") B28=(3016) C28=(2985) B29=(100) C29=(100) ///
A31=("Treatment means are regression adjusted") A32=("*** p<0.01, ** p<0.05, * p<0.1")

//Export variable label, regression coefficient, control mean, treatment mean, p-value, sample size
local outcomes pk07A pk07a_A pk08A pk08aA pk09A pk10A pk11A pk13A pk14A pk15A pk15aA ///
pk07_sample pk07a_sample pk08_sample pk08a_sample pk09_sample pk10_sample pk11_sample pk13_sample pk14_sample pk15_sample pk15a_sample ///
non_dignified_care provider_neglect provider_abuse
local row = 2
foreach var of local outcomes {
	sleep 2000
	local varlabel : var label `var'
	
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	
	local p = 2*ttail(e(df_r), abs(_b[treatment]/_se[treatment]))
	local stars
		if `p' < 0.10 local stars = "*" 
		if `p' < 0.05 local stars = "**" 
		if `p' < 0.01 local stars = "***"
	local impact = _b[treatment]
	
	qui sum `var' if treatment == 0
	putexcel A`row' = ("`varlabel'")
	putexcel B`row' = ((_b[treatment]) + (r(mean)))
	putexcel C`row' = (r(mean))
	putexcel D`row' = ("`impact'" + "`stars'")
	putexcel E`row' = (`p')
	putexcel F`row' = (_b[treatment]/r(sd))
	putexcel G`row' = (e(N))
	local ++row
}

********************************************************************************

* CATEGORY 9. Improved facility cleanliness


***Output Tables***

*Label variables*
lab var pk12A "Recent visit: Cleanliness (1-4 scale)"
lab var pk12_sample "Sample Puskesmas: Cleanliness (1-4 scale)"

*Produce Excel output

//Create the table shell
putexcel set "intermediate_outcomes.xlsx", sheet("Category 9") modify 

putexcel A1=("Outcome") B1=("Treatment Mean") C1=("Control Mean") D1=("Impact") ///
E1=("p-value") F1=("Effect Size") G1=("Sample Size") A5=("Number of Respondents") ///
A6=("Number of villages") B5=(3016) C5=(2985) B6=(100) C6=(100) ///
A8=("Treatment means are regression adjusted") A9=("*** p<0.01, ** p<0.05, * p<0.1")

//Export variable label, regression coefficient, control mean, treatment mean, p-value, sample size
local outcomes pk12A pk12_sample
local row = 2
foreach var of local outcomes {
	sleep 2000
	local varlabel : var label `var'
	
	qui reg `var' treatment strata1 strata2 strata3, vce(cluster iddesa)
	
	local p = 2*ttail(e(df_r), abs(_b[treatment]/_se[treatment]))
	local stars
		if `p' < 0.10 local stars = "*" 
		if `p' < 0.05 local stars = "**" 
		if `p' < 0.01 local stars = "***"
	local impact = _b[treatment]
	
	qui sum `var' if treatment == 0
	putexcel A`row' = ("`varlabel'")
	putexcel B`row' = ((_b[treatment]) + (r(mean)))
	putexcel C`row' = (r(mean))
	putexcel D`row' = ("`impact'" + "`stars'")
	putexcel E`row' = (`p')
	putexcel F`row' = (_b[treatment]/r(sd))
	putexcel G`row' = (e(N))
	local ++row
}

********************************************************************************

clear



